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A phase soliton is a topological defect peculiar to two-band superconductors, which is associated 
with a 2tv winding of the relative phase of the two superconducting condensates. We study the 
quasiparticle spectrum in the presence of a single planar phase soliton. We show that the order 
parameter phase variation in each of the bands leads to the existence of subgap states bound to 
the soliton. Calculation of the soliton energy valid at all temperatures is presented, with exact 
■ analytical results obtained for a simple soliton model. 
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I. INTRODUCTION 



The recent resurgence of interest in the properties of multiband, in particular two-band, superconductors has been 
— — i" largely stimulated by the discovery of superconductivity in MgB 2 (Refs. [l] and0). Other candidates for multiband 
superconductivity include nickel borocarbides (Ref. 0) , NbSe2 (Ref. 0) , the heavy- fermion compounds CeCoIns (Ref. 
H) and CePtaSi (Ref. 0), and also the whole family of iron-based high-temperature superconductors (Ref. 0). These 
discoveries have shown that multiband superconductivity, which is characterized by a significant difference in the order 
parameter magnitudes and/or phases in different bands, is a much more common phenomenon than was previously 
rj \ thought. 

Theoretically, a two-band generalization of the Bardeen-Cooper-Schrieffer (BCS) theory was introduced in Ref. H|. 
Subsequent work has shown that many properties of multiband superconductors differ qualitatively from the single- 
band case and that the most spectacular features are associated with the presence of additional degrees of freedom - 
the relative phases of the pair condensates in different bands. For example, in a charged two-band superconductor, the 
collective mode corresponding to small oscillations of the relative phase, called the Leggett mode,— is not accompanied 
by the charge density modulation and, therefore, is not pushed up into the plasma frequency region. If the two 
condensate phases have different windings around the core of a vortex, then the vortex will carry a fractional magnetic 
O- flux^ 

In addition to exotic vortices, there is another type of topological defects specific to multiband superconductivity, 
, namely phase solitons^i A phase soliton is a topologically stable texture of the superconducting order parameter, in 
which the relative phase exhibits a kink-like variation by 2ir between its asymptotic mean-field values. The phase 
■^j- | solitons can be dynamically generated in nonequilibrium current-carrying states^ or even in static situations by the 
0^ ■ proximity effect with a conventional s-wave superconductor^ Phase solitons of a different kind, connecting degenerate 
time-reversal symmetry breaking states, may exist in superconductors with three or more bandsJ^ Stable nontrivial 
phase textures similar to the phase solitons can also exist in single-band superconductors with unconventional multi- 
component order parameters. For instance, a chiral p-wave superconductor with k x ± ik y gap symmetry can break 
up into domains of opposite chirality separated by a domain wall, in which the relative phase of the order parameter 

■ components rotates between — tt/2 and ir/2 (Ref. [Pot ). 

Previous studies of the phase solitons in multiband superconductors focused on finding the soliton shape and energy 

■ in the Ginzburg-Landau regime* ^ 13 ' 16 In this paper, we investigate the effect of the phase solitons on the Bogoliubov 
quasiparticles. The presence of a nonuniform texture in the relative phase implies that the order parameter phases in 
individual bands also have kink-like inhomogeneities. We show that, in addition to the gapped quasiparticles in the 
bulk, there are states in both bands which are localized near the soliton and have energies below the bulk gap edges. 
The origin of these states is similar to that of the fermion states bound to topological defects, which have appeared 
in many different contexts in high energy and condensed matter physics, 1 - A different type of subgap states that can 
exist near the surface of a two-band superconductor of the s± symmetry, i.e. when the order parameters in the bands 
have opposite signs, was discussed in Ref. [H. 

The paper is organized as follows: In Sec. UH we discuss the structure of the phase soliton in the Ginzburg-Landau 
regime and show how the kink in the relative phase is translated into kinks in the individual condensate phases, with 
non-universal phase winding numbers. In Sec. IIII1 we study the quasiparticle spectrum in the presence of a single 
planar soliton using semiclassical, or Andreev, equations and calculate the energy of the bound states. In Sec. IIV1 
the phase soliton energy is calculated using an exact representation of the functional determinant of the Andreev 
Hamiltonian. Throughout the paper we use the units in which h = ks = 1. 
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II. GINZBURG-LANDAU DESCRIPTION OF THE PHASE SOLITON 



We assume a clean superconductor with two isotropic bands, labeled by a — 1, 2, and isotropic s-wave singlet pairing, 
described by two order parameters rji(r) and ^(r), in the absence of a magnetic field. The difference between the 
free energies in the superconducting and normal states is given by J- s — J- n = j Jgl d 3 r, where 



GL 



E 



a a \fla\ 



Pa I ,4 

yKI 



K a \Vr, a \ 



(1) 



The intraband terms have the usual Ginzburg-Landau form, while the last term describes the interband "Josephson 
coupling" , i.e. the Cooper pair tunneling between the bands. 

Using the amplitude-phase representation of the order parameter, i] a (r) — \i] a {''')\e lLPa '^ r \ the free energy density 
can be written as 



Igl — 



ValVal 2 + y 1^ | 4 + K a (V\ Va \) 2 + K a \ Va \ 2 ( V^ a ) 



' 2-y | r/i 1 1 772 1 cos(ipi - (p 2 ). 
9q (mod 2it) , where 



(2) 



In a uniform state, the minimum energy corresponds to |ry | — A a and (pi — <p2 

O = O, if 7 <0, 9 = it, if 7 > 0. (3) 

The first possibility (interband attraction) is realized in MgB2, in which both gaps have the same phased while the 
second possibility (interband repulsion) is likely realized in the iron pnictides, in which, according to the most popular 
model, the gap function reverses its sign between different sheets of the Fermi surface, corresponding to the so-called 

• • 2D 

s± pairing^ 

It follows from Eq. @ that the supercurrent is a sum of independent contributions from individual bands: j = 
— (4e/c) ^ a K a \r] a \ 2 (Vip a ) (e is the absolute value of electron charge). For a planar texture perpendicular to the x 
axis, the current conservation implies that j = jx, where j is a constant. The value of the current is set by external 
sources and can be assumed to be zero. In order for the supercurrent contributions from bands 1 and 2 to cancel 
each other, the two order parameter phases must vary in a counterphase fashion, with = — p{x)^ xfi, where 

p = Ki\r]i\ 2 / K2\rj2\ 2 ■ This allows one to express the free energy @ in terms of the relative phase 9 = ipi — if2'. 
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(4) 



Variational minimization of this expression yields a system of three coupled nonlinear differential equations for (771 (x)|, 
|?72(a;)|, and 6(x), with the asymptotics |^ a (±oo)| = A a and 9 (±00) = 6* (mod27r). 

In addition to the uniform solutions, the order parameter equations have various nonuniform ones, connecting 
different degenerate minima of cos 9. The simplest topologically nontrivial solutions are those with 9(+oo) — 9(—oo) = 
±2tt, where the positive (negative) sign corresponds to a phase soliton (anti-soliton). The presence of a soliton texture 
in the relative phase implies that each of the two phases ipi and ip2 is also spatially nonuniform and attains different 
values at x = +00 and x = —00. We define the phase winding parameter as 



X = <^i(+oo) - <P\{-oo) 



+00 



dx 



1 + p{x) ' 



(5) 



then (/3 2 (+oo) — <y9 2 ( — 00) = \ T 27r for the soliton (anti-soliton). 

An explicit expression for the phase soliton can be obtained in the London approximation, when the order parameter 
amplitudes are constant everywhere, i.e. |ry (a;)| = A a (Ref. [ill ). The minimization of Eq. (j4]) then yields a static 
sine-Gordon equation for the relative phase, whose soliton solution has the form 9{x) = 9 s (x) + (n — 9q), where 
9 s (x) = 2 arcsin[tanh(a;/£ s )] and 



! K 1 K 2 A 1 A 2 
(K 1 A 2 + K 2 A 2 )\ 1 \ 



has the meaning of the soliton width. The phase textures in the bands are given by the following expressions (up to 
a common phase rotation): 

1 



ipxix) 



Po 



,(x) - (it - 6» ), 



(6) 



1 + Po " ' 1 + Po 

where po — K\ A 2 /K2A 2 . In the London approximation, the phase winding parameter, see Eq. (JSJ), takes the form 
X = 27r/(l + p ). 
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III. QUASIPARTICLE SPECTRUM 

The qualitative features of the phase soliton discussed above are expected to survive beyond the Ginzburg-Landau 
regime. Namely, the phase soliton divides the superconductor into two domains, separated by a "domain wall" , whose 
thickness is of the order of £ s . The order parameter phase in each of the bands exhibits a kink- like variation, similar 
to the London-limit expressions, see Eq. ((6]), with ip a (+oo) — <p a (— oo) — \a- For a single soliton, we have 

Xi = X, X2 = X - 2?r (7) 

where the phase winding parameter x is a non-universal fraction of 2-7T, determined by the microscopic details. 

Now we turn to the calculation of the quasiparticle spectrum in the presence of a single planar soliton. The bands 
are isotropic, with the dispersions £, a {k) = (k 2 — k Fa )/2m a , characterized by the effective masses m a and the Fermi 
wave vectors kp. a - Since the order parameters vary slowly on the atomic length scales, one can use the semiclassical, 
or Andreev, approximation^ An important point is that the slow perturbation due to the phase soliton cannot cause 
quasiparticle transitions between the bands, therefore one can solve the Andreev equations independently in each 
band. 

Quasiparticles propagating along the semiclassical trajectory directed along the unit vector k F are described by 
the wave function -0, which varies slowly compared to kp . The quasiparticle spectrum at given kp is determined by 
the equation Hip = Eip, where the Andreev Hamiltonian is given by 

6=(- iv fr7 x )■ (8) 

y T) (X) lV F ,xVx J 

Here v F = kp/m is the Fermi velocity and r](x) = \ri(x)\e lv ( x \ The gap magnitude approaches its bulk mean- 
field value far from the soliton: \i](x)\ — > A at |x| 3> £ s [the London approximation corresponds to \r](x)\ = A 
everywhere]. While the band index has been temporarily dropped for brevity, we note that in the ath band, v F — > 
VF,a = (kF,a/m a )k F , A -)• A a , and <p(x) -> ip a (x). 

To make the eigenvalue problem for the Hamiltonian © well-defined, we put the system in a box of length I, 
such that I ^> £ s . When safe to do so, we will take the limit t — > oo. For consistency with the phase winding of 
the order parameter, one should use twisted boundary conditions for the quasiparticle wave functions: ip(+£/2) = 
e ix&3 / 2 ip(-e/2). 

It is convenient to represent the phase soliton as a localized perturbation, which is achieved by applying a gauge 
transformation: ip — Uip and U^HU — H, where 

U(x) = e M*)&s/2, (9) 
We can drop the tildas and write the transformed Hamiltonian in the form 

H = H Q + SH, (10) 

where 

H = -iv F<x a 3 V x + A <Ti (11) 
describes the Bogoliubov quasiparticles in the uniform superconducting state, while 

SH = -v FjX cp'(x)a + [\t}(x)\ - A ] u\ 

represents a perturbation which is nonzero only near the soliton, i.e. at |ce| < The gauge-transformed eigenfunctions 
satisfy the periodic boundary conditions: 

*(+{) -*(-£). 

Note that there is a one-to-one correspondence between the spectra of H and Ho: the eigenvalues of the operator 
H s = Ho + ASH evolve smoothly between those of Hq and H as the parameter A varies between and 1. 

At given kp, the spectrum of the Andreev Hamiltonian consists of scattering states with the energies \E\ > Ao and 
bound states with \E\ < A . Let us start with the former. Far from the soliton, the Hamiltonian is equal to H and 
the scattering eigenstates are the superpositions of plane waves: 

± I W R i Jqx i n± ( w l 



^)\^±oo = C R l ^ + G t T e"^, (13) 
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where q = ' E 2 — Aq/|«j? iX | > 0, w_r(l) = ^o/(ET v F, x q), and the subscripts R, L refer to the direction of propagation 
of the corresponding waves. The coefficients in these asymptotics are not independent: it is convenient to introduce 
a 2 x 2 scattering matrix, or the S-matrix, which expresses the amplitudes of the outgoing waves in terms of the 
amplitudes of the incoming waves: 

The elements of the S'-matrix depend on the energy and are determined by the details of the order parameter 
at \ x \ Cs- Note that the S'-matrix defined by Eq. (fT4|) is not unitary, in general, since we did not bother to 
normalize the scattering states. Still, one can show that the S-matrix satisfies a certain constraint, which follows 
from a "conservation law" for the Andreev equations. It is straightforward to check that Vx^a^tp) = f° r the 
eigenfunctions of Eq. (|8]), therefore, ip^ (x)&3tfj(x) = const. Substituting here the asymptotical expressions (1131) and 
using the definition (fl4| . we obtain that the S-matrix must satisfy S^jlS — /t, where fi = diag(wfj — 1, 1 — w\). In 
particular, |detS| = 1. 

One can also introduce the T-matrix, which relates the scattering wave amplitudes at x — >■ +oo to those at x — >• — oo: 

I Cr\ - ( C R \ 



ci \c- L 



(15) 

Comparing Eqs. (|15j) and (|14|) . we find that the r-matrix can be expressed in terms of the S-matrix: 

(16) 

In the absence of the phase soliton, there is no scattering and S = f = (Xq. 




A. Subgap bound states 



The S-matrix (or the r-matrix) can also be used to obtain the bound states, which correspond to the poles at 
\E\ < Aq on the real axis in the complex energy plane. The function q{z) = y/ z 2 — Aq/\vf,x\> where z is the complex 
energy, has two branch points at z = ±Ao. The appropriate branch of q(z) is fixed by the condition that, as implied 
by Eq. (|13[) . q is a positive real number when z is outside the gap on the real axis, i.e. when z = E with \E\ > Aq. 



One can select the branch cuts to run parallel to the imaginary axis, from ±Ao to ±Ao =F ioo. Then, at \E\ < Aq we 
have q — ifl/\vF,x\, where ft — ^/Aq — E 2 . 

The S-matrix can be calculated analytically in a simple model, in which the soliton width is sent to zero, so that 



|r?(a;)|=Ao, <p(x < 0) = 0, <p(x > 0) = x , 



(17) 



where \ is the phase winding parameter. The gauge transformation operator U (x), see Eq. (|9]), is discontinuous at x = 
0, which implies the following matching condition for the gauge-transformed wave function: ip(+0) = e-^ a ^ 2 ip(-0). 
After a straightforward calculation, we obtain: 



/ X , „• E ■ X 
S = I cos — + i sin — 

2 v F , x q 2 



-i 1 



1 

E 

sin -f- 

vf,xQJ 2 



i I 1 sin — 1 

VF.xqJ 2 
l 



(18) 



The characteristic equation for the bound states at \E\ < Aq has the form 



y E \ 

cos - + sgn \v FtX )— sin - = 0. 



(19) 



Introducing E = E sgn (vp,x)i one can write E — AoCOsO and Q = AosinO. Since, according to Eq. (|19j) . 
tan O = — tan(x/2), we have O = — x/2 + irn (n is an integer) and, therefore, E = Aq(— 1)™ cos(x/2). The parity of 
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n can be found from the condition f2 > 0, which yields (—1)™ = — sgn [sin(x/2)]. Collecting everything together, we 
obtain: 

E = -A sgn (v FtX sin ^ cos ^, (20) 

i.e. there is a single bound state with the energy inside the bulk gap. In the absence of the soliton, i.e. at x = 0, 
we have \E\ — Aq, i.e. the bound state merges into the continuum of the bulk states. Note that the "sharp" phase 
soliton is formally similar to a Josephson junction between two s-wave superconductors, with the phase difference 
equal to \, The bound state energy for such a junction was calculated in Ref. I22I . 

Restoring the band indices and using the phase windings from Eq. ((7J|, we finally obtain that there is one subgap 
bound state for each direction of semiclassical propagation kp in each of the bands, with the energy given by 

y 

E a = -A a sgn (v FyayX ) COS ~. (21) 

We see that the bound state energy is a non-universal fraction of the bulk gap. It is only in the exceptional case when 
the microscopic parameters are fine tuned to yield X = 7r, that the subgap states are located exactly at zero energy. 

Expression (|21| has the property E a {— kp) = —E a (kp), which is a consequence of the "electron-hole" symmetry 
of the Andreev spectrum: for any hp, if ifi^ is an eigenfunction of the Andreev Hamiltonian -ff^ F corresponding to 
the eigenvalue E, then ia%ip% is an eigenfunction of the Andreev Hamiltonian H_^ F corresponding to the eigenvalue 
— E. After angular averaging over the Fermi surface, the bound states will manifest themselves as four ^-function 
peaks in the quasiparticle density of states, located symmetrically at E = iAi^ cos(x/2). 



IV. ENERGY OF THE PHASE SOLITON 



Since the quasiparticles bound to the phase soliton have lower energies than in the bulk, it is natural to ask 
whether the spontaneous formation of solitons, accompanied by "self-trapping" of quasiparticles, could be possible. 
The general expression for the energy of a nonuniform state in a two-band superconductor is derived in the Appendix. 
For a planar order parameter texture with rj a {x) — \rj a {x)\e lVa ^ x \ in particular, for the phase soliton, the free energy 
difference per unit area between the states with and without the soliton has the form F s = 5J r /A±, where ST is given 
by Eq. (|A6[) and A± is the area of the system in the directions perpendicular to x. We have F s = Fi + F2, where 



and 



d 2 k± ^ iuj n - E a , k± ^ 



(2tt) 2 ^ 



E 



(0) 
a,h± 



F 2 = J dx^2(V l ) ab (77*77;, - ??* !0 ?7m) 

ab 



(22) 



(23) 



In F\, we used the following notations: uj n = (2n + 1)ttT is the fermionic Matsubara frequency, k± = (k y , k z ) is the 
wave vector parallel to the soliton, and [i labels the eigenstates of the reduced one-dimensional BdG Hamiltonian in 
the ath band, at given fej^: 



frBclG 
n a.k 1 



/ P - k 2 
2m a 



k 2 - h 2 



- k 

2m a ' 



(24) 



In i*2, V a b are the coupling constants of the intraband and interband pairing, see Appendix 



where fco, a ; 
for the explanation. 

Since Fx is a sum of the independent contributions from the two bands, we can drop the band index temporarily. 
As in Sec. IIII1 one can use the Andreev approximation to find the spectrum of the Hamiltonian (IM|) . because the 
order parameter varies slowly on the scale of the inverse Fermi wave vector. We seek the eigenfunctions of Eq. (f24|) in 



the form ^(x) = e * x ip(x), where k x — ±fco- The direction of semiclassical propagation of quasiparticles is defined by 
the wave vector kp = (fej_, k x ) = kpkp. The slowly-varying function ip(x) is found by solving the eigenvalue equation 



6 



Hip — Eip, where H is the Andreev Hamiltonian, see Eq. ((8]). The sum over the BdG spectrum in Eq. (|22j) can be 
expressed in the semiclassical approximation in terms of a Fermi-surface angular average of a sum over the Andreev 
spectrum: 

where Np = mkpjl-K 1 is the Fermi-level density of states and i labels the eigenstates of the Andreev Hamiltonian at 
given ItF. 

Removing the order parameter phase by the gauge transformation (j9]) and restoring the band indices, we finally 
arrive at the following result : 

Ft = -2nTj2J2 NF <-[ ^\ VF >«>*\ ln7 \fc F (^«), (25) 

n a 

where 

tt z - E t (a,k F ) _ dct[z - H(a,k F )] 

U a.k F \ Z ! - 11 _(()), j \ i i\ tt i j M ^ ' 

t z - El ' (a, k F ) det[z - H (a, k F )\ 

is the ratio of the functional determinants of the Andreev Hamiltonians in the nonuniform and uniform states, see 
Eqs. (fT0|) and (fTTj). for a given direction of the semiclassical propagation on the Fermi surface in the ath band. 



A. Calculation of the functional determinant 



There exists a very efficient way of calculating the expression (f2lj)) . which is based on a relation between the 
functional determinant and the transfer matrix for the Andreev Hamiltonian (see Ref. l23l where a closely related 
Dirac Hamiltonian was investigated). Let us again drop the band and direction indices, a and kp. The transfer 
matrix is defined as a 2 x 2 matrix satisfying the equation (z — H)M(x; z) = 0, where H is given by Eq. (fTTH) , with 
the initial condition M(—£/2; z) — ctq. Since the eigenfunctions of H can be written as ip(x) = M(x; z)ip(—£/2), the 
transfer matrix has the meaning of the evolution operator of the wave functions along the x-axis. 

From the periodic boundary condition (|12[) we obtain the characteristic equation for the eigenvalues: det[<To — 
rh(z)] = 0, where we introduced a shorthand notation, m(z) = M(£/2; z), for the transfer matrix from one end of 
the system to the other. The ratio of the functional determinants, Eq. (|2"o]> . can then be represented in the following 
form: 

= y-f(*)] (27) 
det[<r - m (z)\ 

where mo is the transfer matrix from — 1/2 to £/2 for Hq- That the two sides of Eq. (|2T[) have to be the same 
immediately follows from the fact that they both are meromorphic functions in the complex energy plane, having the 
same poles and zeros and also the same asymptotics at \z\ 3> Aq. If this argument is not convincing, a more elaborate 
proof can be found in Ref. HH An expression like Eq. (|2T|) represents a significant step forward compared to the 
definition (|26|) . because it reduces the calculation of the infinitely-dimensional functional determinant to solving an 
initial value problem for a 2 x 2 transfer matrix. Expressions of this sort are sometimes called the Gelfand-Yaglom 
formulas, see Ref. [13 and also Ref. '25] for a review. 

Further simplification is possible in the thermodynamic limit, £ — > oo, where one can represent Eq. (I27|) in terms of 
the scattering matrix. To obtain this representation, we note that, according to the definition of the transfer matrix, 

^ (+^\ = mip ( ~Y (28) 



2) V 2 . 

On the other hand, using Eq. (|13p . the wave functions far from the phase soliton can be expressed in terms of the 
scattering wave amplitudes as follows: 

*&)-**[%)• (2S) 
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where 

W± 



e ±iqi/2 

It follows from Eqs. (fl5 ]> . ([28]> . and ([29]) that m = W+tWZ 1 and 

Diz) = (30) 
det(<r - W+WZ 1 ) 

Here we used the fact that f = &o for Hq. 

According to Eq. (|25|) . the free energy of the phase soliton is expressed in terms of the Andreev functional 



determinant on the imaginary energy axis. At z = iui n , we have q = in, where k = \/uj^-\-A^/\vf,x\- Calculating 
the 2x2 determinants on the right-hand side of Eq. (|3"0")) and keeping only the leading, exponentially divergent at 
I — > oo, terms, we obtain: det(<7 — W+fWI 1 ) — —e Kt T 2 2 and det(<5o — W+iyr 1 ) = —e Kl . Therefore, 

D ( itJ ™)|^oo = T 220n) = TT-^ s , (31) 



where we used the relation (|16p between the r- and S*-matrices. 
Returning to Eq. (|2"5|) , we finally obtain: 

Ft =2TrT^2^N Fia j ^L\v F ^ x \ In S 22 (iu n ; a, k F ). (32) 

n a J 

Thus, the problem of evaluating the free energy of a nonuniform order parameter texture has been reduced to the 
calculation of the semiclassical scattering matrix of the Bogoliubov quasiparticles, analytically continued to complex 
energies. 

B. Sharp soliton 

The scattering matrix can be calculated explicitly only in some simple cases. For instance, for the sharp phase 
soliton defined in Sec. MI A[ it is given by Eq. (fT5)l and we have 



T ln S22 K; a, k F ) = -T J2 In ( 1 ~ 2 ,\ 2 sin 



A a „• 2X 

According to Eq. (|23p , for the sharp soliton F 2 vanishes and we obtain the following exact expression for the energy, 
which is valid at all temperatures: 

F s = n J2 N Fa v F , a T^lnfl- sin 2 f\ . (33) 

a n>0 \ n a / 

At T = 0, the Matsubara sum here becomes an integral and can be calculated in a closed form: 

r duJ 1 (a A ° • 2 X\ A a / X n 

/ —m l 5 T-^T sm 7T = 1 ~ co S 7T 

i 2tt V ^ 2 + A a 2 y 2 V 2 J 

Therefore, 

F s (T = 0) = J (l - cos || ) N F , a v Fta A a . (34) 



Expressions (|3"3")l and (jM]) show that the soliton energy is positive, vanishing only in the absence of the phase 
winding, i.e. at x — 0. Thus we come to the conclusion that the spontaneous formation of the phase solitons is 
energetically unfavorable. Note though that a definitive answer would require a self-consistent solution of the gap 
equations. The feedback effect of the subgap states on the order parameter profile might be strong enough to cause 
self-trapping of Bogoliubov quasiparticles, similar to that discussed in Ref. [2|| see also Ref. H3. Investigation of this 
possibility is beyond the scope of the present work. 
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V. CONCLUSIONS 



We studied the Bogoliubov quasiparticle spectrum in a two-band superconductor, in the presence of a soliton- 
like topological defect in the relative phase ipi — (p 2 . While the relative phase winding across the soliton is given 
by 2-7T, the phase windings in individual bands are non-universal fractions of 2tt: </? 1 (+oo) — ifi (—00) = x an d 
</?2(+oo) — V>2{— 00) = X — 27r, where the parameter \ depends on the microscopic details. We found that there are 
quasiparticle bound states localized near the soliton, whose energies are non-universal fractions of the bulk gaps. 

The bound states will lead to sharp peaks in the quasiparticle density of states at E — iAj.,2 cos(x/2), which can be 
observed in tunneling experiments. The tunneling probe will have to be located sufficiently close to the phase soliton 
to be able to detect the contribution from the localized states. This can be done, e.g. in the experimental setup 
proposed in Ref. [l3l in which the soliton is "pinned" to the spatial variation of the interband Josephson coupling, 
controlled by the proximity effect with another superconductor. We note that the peaks in the density of states are 
expected to acquire a finite width when impurity scattering or intraband gap anisotropy are taken into account. 

We also derived a general expression for the phase soliton energy, relating it to the scattering matrix of the 
Bogoliubov quasiparticles. As a simple application, we exactly calculated the energy in the limit of zero soliton width. 
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Appendix A: Free energy of a nonuniform two-band superconductor 

In this Appendix, we present a microscopic derivation of the free energy of a clean two-band superconductor 
in a nonuniform state, at arbitrary temperature, using the effective action formalism. We start with a two-band 
generalization of the BCS Hamiltonian: 

# = £fo(*04,a,oA.o,a- y Yj V ^i+q. A .^-k.a.,i d -k' Ai^k'+qM^ (Al) 
k,a kk'q,ab 

where a — 1, 2 is the band index, a =t, \. is the spin projection (the spin indices that appear twice are summed over), 
and V is the system volume. The second term in the Hamiltonian describes singlet s-wave pairing interactions in 
the Cooper channel: the intraband pairing, characterized by the coupling constants V\\ and V22, and the interband 
"tunneling" of the pairs, described by V12 and V21. The Hermiticity and time-reversal invariance of the Hamiltonian 
dictate that the constants V a b form a real symmetric matrix V . 

Our derivation is a straightforward generalization of the standard textbook procedure in the single-band case, see, 
e.g., Ref. The partition function for the Hamiltonian (| Al|) can be represented as a functional integral over the 
Grassmann fields Cfc. a . Q (r) and Cfc iaiQ (r): Z = J DcDc e~ s ^ c ' c \ where the action is given by 



rP Q rP 

S= drS^ c k . a . a — c kaa + / drHl^c], 
Jo tl dT Jo 



with /3 = 1/T. The interaction term in the action can be written as 



S int = - f dTVy2V ab B a (q.T)B b (q.T), 

JO u 
q.ab 



where 

Ba{q, T) = y ^ Cfe+q,a,t( T )C-fc,a,|(r), B a (q, t) = y ^ C-k^i^Ck+g^T). 
k 



One can decouple the interaction by means of the Hubbard-Stratonovich transformation, introducing two complex 
conjugated bosonic fields 771,2(9, t): 



e -S in t _^ 



Vr,*V V exp J - f dr^J^iV-'Uvlrjb- f dr YfaB* + S aVa) \ 

[ J ° V q,ab J ° q,a J 
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The field rj a has the meaning of the fluctuating order parameter in the ath band. 

One can now calculate the Gaussian integral over the fermionic fields, to obtain Z = J Vrj*Vri e~ Se/ * f [' ) ' v \ where 



Vb (A2) 



is the effective bosonic action and 



G~ = 



-dr-^a(k) -Va(r,r) 
-r?*(r,T) -d T + Uk) 



is the inverse Green's operator in the ath band, with k = — iV. The trace in the first term should be understood 
as an operator trace in (rr)-space and a 2 x 2 matrix trace with respect to the electron-hole (Nambu) indices in the 
ath band. Near the critical temperature, the order parameter components are small and the first term in Eq. ()A2|) 
can be expanded in powers of r\ a . In this way one would arrive at the Ginzburg-Landau functional for the two-band 
superconductor, which has been derived by different means in Ref. [29l 

We do not restrict ourselves to the Ginzburg-Landau regime and calculate the free energy in the mean-field approxi- 
mation at arbitrary temperature. The mean-field solution for the order parameter corresponds to a static saddle point 
of the effective action and satisfies the equations SS/Sr)* 2 = 0. From Eq. (|A2j) we obtain two coupled self-consistency 
equations: 

%(»•) = -T^J]v; 6 G M2 (r,r;w„). (A3) 

n b 

Here uj n = (2n + l)nT is the fermionic Matsubara frequency and G a ^i is the anomalous (Gor'kov) component of the 
matrix Green's function G a , which satisfies the equation (iu) n — H a 3dG )G a (r,r';uj n ) = &oS(r — r'), where 

frBdG _ ( ta{k) ri a {r) \ ( 
Ha -U(r)-6,(*)J (A4) 

is the Bogoliubov-de Gennes (BdG) Hamiltonian. 

The gap equations can be represented in terms of the eigenstates and eigenvalues of the BdG Hamiltonian, which 
are found from H^ dG ^ ap (r) — E a:P ^ a ^ p (r). The eigenstates are two-component Nambu spinors, ^ = (u, v) T , labeled 
in the ath band by quantum numbers p. Assuming that the eigenstates form a complete and orthonormal set, we 
obtain for the matrix Green's function: 



G a (r,r';uj n ) = 



j2 ^-Ar)HA r ') 

p 



Inserting this into Eq. (|A3I) and using the "electron-hole" symmetry of the BdG spectrum (if \& corresponds to the 
energy E, then z^^* corresponds to the energy — E), we arrive at the final form of the gap equations: 

Va(r) =E^E ' u bArK, p (r)[l - 2f(E b>p )}, (A5) 

b p 

where f(E) = l/(e^ E + 1) is the Fermi function. The prime means that the summation is performed only over the 
upper half of the BdG spectrum, i.e. over the eigenstates with Eb :P > 0. Eq. (|A5[) can be used to obtain both the 
critical temperature and the temperature dependence of the gaps, see Refs. H and In addition to the spatially 
uniform solution, given by riifi — Aie* e °, 772,0 = A2, where 9q — for interband attraction (V12 > 0) and 6q = n 
for interband repulsion (V12 < 0), the gap equations also have various nonuniform solutions, in particular, the one 
corresponding to the phase soliton. 

The effective action (|A2I) for any mean-field configuration of the order parameter has the form S e ff = f3£, where 

£ = -7^]Tln( lW „-£ a , p )+ I d'rJ^iV-XbVlMvbir). 

n a,p ab 



The mean-field free energy is given by T — —TlnZ — const + £. To remove the undetermined constant, we calculate 
the free energy difference between a given nonuniform superconducting state and some reference state. For our 
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purposes, it is natural to choose the latter to be a uniform superconducting state with the order parameters equal to 
r] a fi, and we finally obtain: 

st = m Hm] = -rEE ln : E 7) + / ^r^iv-'uivfa-nlonw). (A6) 

n a,p ^a,p J a b 

Here E$l are the eigenvalues of the BdG Hamiltonian (|A4[) in the uniform state. 
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